perm filename DHSYIQ.SAI[PIC,HE] blob
sn#430334 filedate 1979-04-03 generic text, type C, neo UTF8
COMMENT ⊗ VALID 00012 PAGES
C REC PAGE DESCRIPTION
C00001 00001
C00003 00002 ENTRY "dhsyiq"
C00008 00003 ! This is a short description of what the algorithms really do.
C00012 00004 ! this page contains all declarations and global procedures
C00015 00005 ! this page contains the macros necessary for the hue calculations
C00018 00006 ! this page contains the code to initialize files and related variables
C00021 00007 ! this page contains the code to initialize various constants
C00022 00008 ! this page contains the code to preload various arrays
C00024 00009 ! this page contains the code for beginning the main outer (rows) and inner (columns) loops
C00026 00010 ! this page contains the guts of the program
C00029 00011 ! this page clears buffers and releases them
C00030 00012 REQUIRE UNSTACK!DELIMITERS
C00031 ENDMK
C⊗;
ENTRY "dhsyiq";
BEGIN "dhsyiq"
REQUIRE "35A" COMPILER!SWITCHES;
REQUIRE "vislib.dcl" SOURCE!FILE;
REQUIRE "⊂⊃<>" DELIMITERS;
! Rewritten by John R. Kender, July, 1975.
Modified by John R. Kender, November, 1976.
This program is a combination and modification of some previous
programs whose purposes were:
1) Transform the original red, green, blue information data into the
psychologically descriptive attributes of intensity, hue, and
saturation.
2) Transform the original red, green, blue information data into the
linear transformations that the color television industry uses.
These are called "Y", "I", and "Q", and correspond roughly to
luminance, orange-red filtration, and magenta filtration,
respectively.
The procedure assumes that the buffer routines have been initialized
by a "bufinit" call outside itself.
Intensity (density) is simply arithmetic mean of the tristimulus
values, rounded.
Hue, saturation and the Y, I, and Q transformations are calculated as
described in the CMU Technical Report, "Saturation, Hue, and Norm-
alized Color: Calculation, Digitization Effects, and Use" by John R.
Kender.
To save time, many of the repeated calculations are done before-
hand, and their values stored in arrays. Thus the intensity, satur-
ation, Y, I, and Q values are calculated using table lookups rather
than by straight multiplication or division. Similarly for part of
the hue calculation.
For ease of programming (mine), the hue calculation is always
performmed assuming the case of a 360 degree output. If a smaller
byte size is necessary, the values are divided down by 3, 5, or 15.
Note to reader: feel free to change this.
Restrictions: The program is designed to handle input files
containing up to and including eight bits per byte. If it is
necessary to process files with a larger byte size, change the macro
"!biggest!byte!handled!" to the new byte size and recompile. Some
reprogramming may be necessary for hue at these byte sizes. Note
also that most of these transformations make little sense if the
input byte size is very small.
Future note: The performance can be greatly improved if some sort of
anticipatory buffering is included (or else, buffering in units of
blocks). Note also that the standard picbuf practice of creating a
file of zeros, then reading in a row, modifying it and rewriting it
is unnecessary. Since the other six files are being created, there
is no need to "update" the new all zero files, and they can be
written directly (say, with a SAIL "arryout" statement).
;
! This is a short description of what the algorithms really do.
Intensity:
i ← (r + g + b) / 3.0 + .5, that is, rounded mean.
As r, g, b integers, integer division with proper prefudge gives
identical result:
i ← (sum + 1) % 3, which is easily converted to table lookup.
Hue:
After deciding which of the six parts of the color triangle to use,
the arctangent is derived as follows:
a ← max - min
b ← mid - min
That is, a and b are saturated color components.
ratio ← (a - b) / (a + b)
This gives relative distance from midpoint of triangle side:
that is, if ratio=0, angle is 0. ratio=1, angle is 60 degrees.
Then, ratio ← sqrt(3) * ratio
gives a proper tangent for the angle under consideration.
So, angle = arctan (sqrt (3) * (max - mid) / (max - min + mid - min)).
The value of the tangent is used as an index into a predictor array:
probe ← atanpredictor[ratio * xfactor]
where xfactor is determined by the size of the predictor array.
tantable[probe] is the tangent of the predicted (integer) angle.
If it is too small, then probe + 1 is used.
Note that xfactor is a constant, hence ratio * xfactor is really
(a - b) * (xfactor / (a + b)). The factor xfactor / (a + b)
can be converted to table lookup.
Note also that the multiplication by the coefficient sqrt(3) can
also be embedded into the tables (as is a factor of 1024
in the case of tantable).
Saturation:
s ← (1 - 3 * (r min b min g) / (r + b + g) )
This gives a saturation in the range 0 to 1. Then,
s ← expander * s + .5
That is, it is plumped up to an integer, rounded.
Rearranging terms:
s ← (expander + .5) - ( (3 * expander) / sum) * (r min b min g)
The term (3 * expander) / sum can be converted to table lookup.
Y, I, and Q:
The matrix multiplication is:
y .509 1.000 .194 r
i = 1.000 - .460 -.540 g
q .403 -1.000 .597 b
As r, g, b integers, each individual multiplication can be
accomplished by a table look up.
Rounding and the proper offsetting of negative values can
be piggybacked into the array lookup, also.;
! this page contains all declarations and global procedures;
INTERNAL PROCEDURE dhsyiq (STRING ARRAY filenames;
BOOLEAN randomize);
BEGIN "dhsyiq!code"
! reassign the following, and recompile, if files are > 8 bits/byte;
DEFINE !biggest!byte!handled! ← ⊂ 8 ⊃ ;
SHORT BOOLEAN continue, no, yes;
SHORT BOOLEAN SAFE ARRAY want[4 : 9];
SHORT INTEGER bytz, colz, f, flag, ii, isrt, jj, jsrt, rowz;
SHORT INTEGER SAFE ARRAY buf[1 : 9], outbytsz[4 : 9];
SHORT INTEGER rp, gp, bp, dp, hp, sp, yp, ip, qp,
r , g , b , d , h , s , y , i , q ;
SHORT INTEGER achromatic, black, satfact, satrange, sum;
SHORT INTEGER hue!divider, hue!range, hue!rounder, hue!achromatic;
SHORT INTEGER cyan, magenta, yellow;
SHORT INTEGER bytemax, bytemaxx2, bytemaxx3, probe;
SHORT REAL angle, deg!per!rad, offset!iq, pi, rad!per!deg,
ratio, sattop, sqrt3;
SHORT REAL rrand, grand, brand, sumrand;
STRING dev;
! declarations used for the various arrays;
DEFINE !biggest!value!handled! ← ⊂ 2 ↑ !biggest!byte!handled! - 1 ⊃ ;
SHORT INTEGER SAFE ARRAY intens [0 : 3 * !biggest!value!handled!];
SHORT INTEGER SAFE ARRAY atanpredictor[0 : 1024 ];
SHORT REAL SAFE ARRAY t1024over [0 : 2 * !biggest!value!handled!];
SHORT REAL SAFE ARRAY tantable [0 : 60 ];
SHORT REAL SAFE ARRAY satfactover [0 : 3 * !biggest!value!handled!];
SHORT REAL SAFE ARRAY yrfact, ygfact, ybfact,
irfact, igfact, ibfact,
qrfact, qgfact, qbfact
[0 : !biggest!value!handled!];
SIMPLE SHORT REAL PROCEDURE tan (SHORT REAL x);
RETURN (SIN (x) / COS (x));
! this page contains the macros necessary for the hue calculations;
DEFINE !arctan! (max, mid, min) ← ⊂ (
IF tantable[
(probe ← atanpredictor[
ratio ← (max - mid)
* t1024over[max - min + mid - min]
]
)
] ≥ ratio THEN
probe
ELSE
probe + 1
) ⊃ ;
! macro for calculating randomized arctan using tables;
DEFINE !arctan!rand! (max, mid, min) ← ⊂ (
IF tantable[
(probe ← atanpredictor[
ratio ← (max - mid)
* 1024 / (max - min + mid - min)
]
)
] ≥ ratio THEN
probe
ELSE
probe + 1
) ⊃ ;
! macro for exploiting the sixfold hue symmetry;
DEFINE !hue!muncher! (red, green, blue, function) = ⊂ (
IF green > blue THEN
IF red > green THEN
yellow - function (red, green, blue)
ELSE
IF red > blue THEN
yellow + function (green, red, blue)
ELSE
cyan - function (green, blue, red)
ELSE
IF green > red THEN
cyan + function (blue, green, red)
ELSE
IF blue > red THEN
magenta - function (blue, red, green)
ELSE
IF blue > green THEN
IF (h ← magenta + function (red, blue, green)) < 360 THEN
h
ELSE
0
ELSE
IF red > green THEN
0
ELSE
achromatic
) ⊃ ;
! macros for dividing the calculated hue down;
DEFINE !divided!hue! = ⊂ (
IF (h ← !hue!muncher! (r, g, b, !arctan!)) = achromatic THEN
hue!achromatic
ELSE
((h + hue!rounder) % hue!divider) mod hue!range
) ⊃ ;
DEFINE !divided!hue!rand! = ⊂ (
IF (h ← !hue!muncher! (rrand, grand, brand, !arctan!rand!)) = achromatic THEN
hue!achromatic
ELSE
((h + hue!rounder) % hue!divider) mod hue!range
) ⊃ ;
! this page contains the code to initialize files and related variables;
! read in red, green, and blue files;
FOR f ← 1 thru 3 DO
BEGIN
dev ← getdev (filenames[f], "dat");
indmp (dev, filenames[f], buf[f] ← fndbuf (0), flag);
END;
! set up switches for files desired;
yes ← -1;
FOR f ← 4 thru 9 DO
IF LENGTH (filenames[f]) ≠ 0 THEN
want[f] ← yes;
! set up other file sizes and get buffers for them;
! note that the red file determines the parameters of the others;
rowz ← rows (buf[1]);
colz ← colms (buf[1]);
bytz ← bytsz (buf[1]);
isrt ← isubst (buf[1]);
jsrt ← jsubst (buf[1]);
! are we properly equipped and outfitted;
no ← 0;
IF bytz > !biggest!byte!handled! THEN
BEGIN
bprmpt ("You have bitten off more than I can chew. See comment in source. Continue? ", continue ← no);
IF NOT continue THEN RETURN;
END;
outbytsz[4] ← bytz;
outbytsz[6] ← bytz - 1;
outbytsz[7] ← bytz + 1;
outbytsz[8] ← bytz + 1;
outbytsz[9] ← bytz + 1;
! find the hue byte size and other related variables;
pi ← 3.1415926536;
FOR hue!divider ← 1, 3, 5, 15 DO
BEGIN
hue!range ← 360 % hue!divider;
hue!rounder ← hue!divider % 2;
outbytsz[5] ← LOG (hue!range) / LOG (2) + 1;
IF 2 ↑ bytz * 2 * pi ≥ hue!range THEN
DONE;
END;
hue!achromatic ← hue!range + 1; ! note: convention only;
! initialize output files;
FOR f ← 4 thru 9 DO
IF want[f] THEN
BEGIN
getbuf (rowz, colz, outbytsz[f], buf[f] ← fndbuf (0));
putsub (isrt, jsrt, buf[f]);
END;
! this page contains the code to initialize various constants;
black ← 0; ! note: convention only (and a bad one at that);
bytemax ← 2 ↑ bytz - 1;
bytemaxx2 ← bytemax * 2;
bytemaxx3 ← bytemax * 3;
deg!per!rad ← 180 / pi;
offset!iq ← 2 ↑ (outbytsz[8] - 1);
rad!per!deg ← pi / 180;
sqrt3 ← sqrt (3);
satrange ← 2 ↑ outbytsz[6] - 1;
satfact ← satrange * 3;
sattop ← satrange + .5;
! code for special hue stuff;
cyan ← 180;
magenta ← 300;
yellow ← 60;
achromatic ← 361;
! this page contains the code to preload various arrays;
! intensity;
FOR ii ← 0 thru bytemaxx3 DO
intens[ii] ← (ii + 1) % 3; ! includes rounding;
! hue;
FOR ii ← 0 thru 60 DO
tantable[ii] ← 1024 * (tan ((ii + .5) * rad!per!deg) / sqrt3);
FOR ii ← 0 thru 1024 DO
atanpredictor[ii] ← deg!per!rad * atan (ii * sqrt3 / 1024) + .5;
t1024over[0] ← 0; ! never accessed;
FOR ii ← 1 thru bytemaxx2 DO
t1024over[ii] ← 1024 / ii;
! saturation;
satfactover[0] ← 0; ! never accessed, as sum = 0 is the singularity;
satfactover[1] ← 0; ! if sum = 1, then (r min g min b) = 0;
satfactover[2] ← 0; ! if sum = 2, then (r min g min b) = 0;
FOR ii ← 3 thru bytemaxx3 DO
satfactover[ii] ← satfact / ii;
! y, i, and q;
FOR ii ← 0 thru bytemax DO
BEGIN
yrfact[ii] ← .509 * ii ;
ygfact[ii] ← 1.000 * ii ;
ybfact[ii] ← .194 * ii + .5 ; ! includes rounding;
irfact[ii] ← 1.000 * ii ;
igfact[ii] ← - .460 * ii ;
ibfact[ii] ← - .540 * ii + .5 + offset!iq; ! includes rounding;
qrfact[ii] ← .403 * ii ;
qgfact[ii] ← -1.000 * ii ;
qbfact[ii] ← .597 * ii + .5 + offset!iq; ! includes rounding;
END;
! this page contains the code for beginning the main outer (rows) and inner (columns) loops;
FOR ii ← 1 thru rowz DO
BEGIN "ii"
rp ← inptr (ii, 1, buf[1]);
gp ← inptr (ii, 1, buf[2]);
bp ← inptr (ii, 1, buf[3]);
IF want[4] THEN
dp ← outptr (ii, 1, buf[4]);
IF want[5] THEN
hp ← outptr (ii, 1, buf[5]);
IF want[6] THEN
sp ← outptr (ii, 1, buf[6]);
IF want[7] THEN
yp ← outptr (ii, 1, buf[7]);
IF want[8] THEN
ip ← outptr (ii, 1, buf[8]);
IF want[9] THEN
qp ← outptr (ii, 1, buf[9]);
FOR jj ← 1 thru colz DO
BEGIN "jj"
r ← ILDB (rp);
g ← ILDB (gp);
b ← ILDB (bp);
IF want[4] OR (want[6] AND NOT randomize) THEN
sum ← r + g + b;
IF randomize THEN
BEGIN
rrand ← r + RAN (0);
grand ← g + RAN (0);
brand ← b + RAN (0);
IF want[6] THEN
sumrand ← rrand + grand + brand;
END;
! this page contains the guts of the program;
! intensity;
IF want[4] THEN
IDPB (intens[sum]
, dp);
! hue;
IF want[5] THEN
IF NOT randomize THEN
IF hue!divider = 1 THEN
IDPB (!hue!muncher! (r, g, b, !arctan!)
, hp)
ELSE
IDPB (!divided!hue!
, hp)
ELSE
IF hue!divider = 1 THEN
IDPB (!hue!muncher! (rrand, grand, brand, !arctan!rand!)
, hp)
ELSE
IDPB (!divided!hue!rand!
, hp);
! saturation;
IF want[6] THEN
IF NOT randomize THEN
IDPB (IF sum = 0 THEN
black
ELSE
sattop - (r MIN g MIN b) * satfactover[sum]
! note that the type of the above expression is
integer, due to the "THEN" part;
, sp)
ELSE
IDPB (IF sumrand = 0 THEN
black
ELSE
satrange * (1 - 3 * (rrand MIN grand MIN brand) / sumrand) + .5
, sp);
! y;
IF want[7] THEN
IDPB (y ← yrfact[r] + ygfact[g] + ybfact[b]
, yp);
! i;
IF want[8] THEN
IDPB (i ← irfact[r] + igfact[g] + ibfact[b]
, ip);
! q;
IF want[9] THEN
IDPB (q ← qrfact[r] + qgfact[g] + qbfact[b]
, qp);
END "jj";
IF (ii MOD 25) = 0 THEN
outst (CVS (ii) & " Rows happily munched so far!" & crlf);
END "ii";
! this page clears buffers and releases them;
FOR f ← 1 thru 3 DO
frebuf (buf[f]);
FOR f ← 4 thru 9 DO
IF want[f] THEN
BEGIN
dev ← getdev (filenames[f], "dat");
outdmp (dev, filenames[f], buf[f], flag);
frebuf (buf[f]);
END;
END "dhsyiq!code";
REQUIRE UNSTACK!DELIMITERS;
END "dhsyiq";